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Waves propagating through a bounded plasma can rearrange the densities of states in the six-dimensional 
velocity-configuration phase space. Depending on the rearrangement, the wave energy can either increase or 
decrease, with the difference taken up by the total plasma energy. In the case where the rearrangement is 
diffusive, only certain plasma states can be reached. It turns out that the set of reachable states through 
such diffusive rearrangements has been described in very different contexts. Building upon those descriptions, 
and making use of the fact that the plasma energy is a linear functional of the state densities, the maximal 
extractable energy under diffusive rearrangement can then be addressed through linear programming. 


I. INTRODUCTION 

Waves propagating through a bounded plasma can 
rearrange the densities of states in the six-dimensional 
velocity-configuration phase space. When the rearrange¬ 
ment is such as to cause particles to diffuse from higher 
energy states to lower energy states, the waves extract 
energy from the plasma. A particular case of this posing 
of the rearrangement problem is the case of alpha chan¬ 
neling, where the energy is deliberately extracted from 
the population of a-particles that are produced in a fu¬ 
sion reactor! 1 This energy is recovered as wave energy. In 
a reactor, this energy is more useful in the form of wave 
energy, which can be used to attain a hot-ion mode or to 
drive electrical current. 

The rearrangements contemplated in plasma using 
waves are diffusive in nature, since the wave-particle 
mechanisms generally cannot maintain coherence, at 
least not for the leading way of using external rf sources 
to heat or drive current in plasma. However, depending 
on the wave frequency and wavenumber, the diffusion oc¬ 
curs in paths that link energy to the spatial dimensions. 
In the case of tokamak reactors, where <a-particles are 
expected to be born at high energy in the plasma center, 
there is a natural energy inversion along the path that 
connects the dense phase space location at high energy in 
the center to the under-dense phase space location that 
is at low energy on the periphery. It is then only a mat¬ 
ter of constructing the appropriate wave diffusion path 
to link these locations. How much energy can be released 
from a-particles is a matter of considerable practical in¬ 
terest, since if appreciable energy could be released in this 
manner, there then might be the opportunity to dimin¬ 
ish substantially the cost of electricity through tokamak 
fusion. 

A theoretical issue of academic interest, however, is 
the precise maximum available energy under diffusive re¬ 
arrangements of phase space density when an arbitrary 
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number of diffusion paths can be constructed. 2 ^ The en¬ 
ergy extractable under diffusion was posed recently, in 
fact, as one of the interesting, outstanding problems in 
wave-particle physics in plasmaP This issue also moti¬ 
vated to some extent approaches to other bounds on en¬ 
ergy exchange between light and plasma, such as the ex¬ 
tent to which bounds could be placed on the absorption 
of laser light at an interface. 4 There are, to be sure, also 
other formulations of free energy in plasma under phase 
space density rearrangements. For example, respecting 
phase space conservation, the Gardner restacking 5 6 rep¬ 
resents a precisely definable free energy that can be read¬ 
ily calculated. This energy is an upper bound on that 
which can be extracted through diffusive processes. How¬ 
ever, because the Gardner restacking is without the real¬ 
istic limitation of the diffusion constraint, the free energy 
available under this formulation represents a rather rar¬ 
efied theoretical construct, even further from practical 
considerations than the academic issue posed here. 

The free energy under the constraint of diffusive re¬ 
arrangements has one well-known textbook example, the 
famous so-called “bump-on-tail” problem. The tail of the 
Maxwellian distribution is imagined to have a “bump” in 
velocity space, so that the distribution is no longer mono- 
tonically decreasing in energy. In that case, waves can 
diffuse particles so as to smooth out the bump, releasing 
the kinetic energy, until a distribution monotonically de¬ 
creasing in energy is reached. The maximum extractable 
energy is obvious, and can be constructed from geometri¬ 
cal considerations; for the 1-bump problem, it is just the 
energy change in flattening the bump. However, were 
there two bumps in the velocity space, then the optimal 
solution would no longer be obvious at all, since it is not 
clear which bump should be flattened first. 

More generally, what is imagined here is that diffu¬ 
sion paths can be constructed that link any two phase 
space locations in the 6D velocity-configuration space, 
whether or not the locations are contiguous.^ The con¬ 
tiguous constraint can formally be realized, in any event, 
in the limit of vanishingly thin paths. Thus we imag¬ 
ine an ensemble of discrete phase space locations, each 
with an initial density, and each representing a certain 
energy. Then as population densities relax under diffu- 
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sion, the total system energy relaxes as well. The free 
energy under the diffusion constraint is then defined as 
the maximum extractable energy, given the opportunity 
to diffuse particles between any two phase space loca¬ 
tions, with any sequence of such two-location or what we 
might call two-state relaxations. It turns out, however, 
that while this pair-wise relaxation is a well-defined pos¬ 
ing of the free energy, it has not been apparent at all how 
to calculate it efficiently when there are many states. 

The problem thus posed in plasmas can similarly be 
posed with respect to stimulated emission by a set of 
lasersP Suppose, an atomic system with just three en¬ 
ergy levels, the ground state at energy ei, the first ex¬ 
cited state at 62 , and the second excited state at 63 , with 
initial population densities of, respectively, nj, rig, and 
Uq . The total energy can be put as IT = e-n, where e and 
n represent the energy levels and the population densi¬ 
ties. Suppose further the availability of three lasers with 
frequencies 1 / 10 , z/ 20 ? and ^ 21 , that, respectively, can stim¬ 
ulate transitions between the first level and the ground 
state, the second level and the ground state, and the 
second level and the first level. Suppose that these lasers 
are incoherent, so what they can accomplish is to equalize 
the populations in any two levels. The maximum energy 
is extracted when the correct sequence of laser pulses is 
applied. 

To make the issues here clear, consider the following 
example. 2 Suppose that the accessible states have ener¬ 
gies with numerical values (0,1,4) and the initial state 
densities are (0,2/7,5/7), where the sum has been nor¬ 
malized to 1. Then the initial energy is Wo = 22/7. The 
energy at step j can be put as Wj. The sequence of 
level-equalizing steps, (^ 21 ,^ 20 ,^ 10 ), then gives 
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The energy extracted is thus 22/7—11/8 = 99/56, or ap¬ 
proximately 56% of the initial plasma energy. One can 
show, for this set of energy levels and initial populations, 
that the sequence used is the optimal sequence for ex¬ 
tracting energy, resulting in the maximum extractable 
energy. What is of interest, however, is how exactly this 
can be proved, how the maximum extraction can be cal¬ 
culated efficiently, and how the complexity of the problem 
increases with the number of states. 

It turns out that the answer to these questions lies in 
the mathematical developments in other fields (although 
these developments appear not to have received much 
attention). Similar level-mixing operations have been 
considered in chemical reaction kinetics.® More directly 
of use here, Zylka identified the set of accessible states 
through level-mixing operations, called the K set, using 
as an example the problem of attainable temperatures 
in heat reservoirs pairwise connected by heat pipes, with 


an arbitrary number of reservoirsP Thon and Wallace, 
in the context of characterizations of altruism as a pair¬ 
wise relaxation correction to economic inequality, derived 
further features of the K setP 

Using these characterizations of the K set, the se¬ 
quence of operations to extract maximal energy can be 
found. Although K, the set of states that can be reached 
from these level-mixing operations, is not convex,® one 
can imagine covering the entirety of K with a small con¬ 
vex polygon. This polygon (namely, the convex hull of 
K, denoted ch(K) here) is determined by taking all pos¬ 
sible convex combinations of the points found within K. 
Equipped now with a linear objective function (the to¬ 
tal energy of a state) and a convex feasible region (the 
unique covering polygon, ch(K )), we may apply the fun¬ 
damental theorem of linear programming to locate the 
minimum energy state at a vertex or edge joining two 
or more vertices of the convex hull of K. Crucially, K 
contains each of the vertices of its convex hull (by con¬ 
struction), and so the minimum energy state over the 
polygon covering K is identical to the minimum energy 
state found within K. 

Thus, using the results of Zylka and of Thon and Wal¬ 
lace, for the three-level system, we shall pose and answer 
the following five questions: 

1. What sequence minimizes the energy? 

2. What is the full set of sequences that must be con¬ 
sidered before the optimal sequence can be found? 

3. What are the full set of sequences that could pos¬ 
sibly be a solution, for some values of the energy 
levels and the initial population densities? 

4. If it were possible to partially relax the distribution 
between two states, rather than fully relax it, would 
that ever be a useful step? 

5. Is it the case that it is ever useful to take a step 
that increases the energy rather than decreases it? 

The paper is organized as follows. In Sec. II, we de¬ 
fine the diffusion model. In Sec. Ill, we reproduce re¬ 
sults of the space of relaxation solutions, and show how 
this immediately answers the first three questions. In 
Sec. IV, we answer the fourth question, proving that par¬ 
tial relaxation is never a useful step, regardless of the 
number of states. In Sec. V, we demonstrate out that 
strategies previously considered 2 ! can be put more pre¬ 
cisely and further prove for any number of states that 
energy-increasing steps can never be part of the optimal 
sequence. In Sec. VI, we offer further discussion of the 
implications for N states and the degree of complexity 
of the problem. In Sec. VII, we summarize the main 
conclusions. 
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II. DIFFUSION MODEL 

Suppose a set of level energies e = (ei, 62 , ..., €n) 
and a set of initial populations no = (nj, rig, ..., riff). 
Here we understand the symbols e and no, without su¬ 
perscripts, to represent vectors; the initial level density 
of state i is represented as n l 0 . Consider a diffusion op¬ 
eration that equalizes the populations of a pair of levels 
(i, j ) such that 



leaving all other level populations unchanged. We are 
interested here is the minimization of the system energy 
Wd = e • rif after repeated application of operations of 
this type on an initial state no and reaching a state rtf. 

The diffusion operation ([T]) can be represented by dou¬ 
bly stochastic matrices of the form = \{I + Qij), 
where / is the TV x TV identity matrix and Qij is the per¬ 
mutation matrix that exchanges the i th and j th level pop¬ 
ulations. Application of equalizes the i th and j th pop¬ 
ulations. For example, 0 could be styled no 
The B^ are symmetric, idempotent, and do not generally 
commute. (Two B^ commute if neither of them operates 
on the same level.) They are a particular case of the T- 
transform: T = (1 — a) I + aQ, a G [0,1] and Q is a per¬ 
mutation of the identity matrix which exchanges only two 
rows.® Like all doubly stochastic matrices, T-transforms 
are measure-preserving: ^2 ab T ab n b = ^ b n b = 1. It is 
important to distinguish between the cases of 0 < a < 
1/2 and 1/2 < a < 1; the latter case corresponds to 
moving density from a less-populated level to a more- 
populated level. For now, we will consider only a = 1/2 
transforms, i.e. the B^. 

III. THE SET OF POINTS K 

Zylka was the first to identify the set of accessible 
states through level-mixing operations: the K setP K 
is determined by the provided set of level populations, 
represented as a length-TV vector no, as well as the al¬ 
lowed diffusion operations. Without loss of generality, 
no may be ordered increasing, as we assume throughout 
this work. As described previously, the diffusive Bij op¬ 
erations equalize the populations of any pair of levels i 
and j. Each element contained in K is an TV-tuple of 
level populations that can be reached by applying some 
sequence of the various B\j to the initial state no- 

The linear function Wd assumes its extremal values on 
the boundary of ch(K). The three-level system is conve¬ 
niently depicted in n\-ri 2 space. Due to normalization of 
the population vector, one coordinate is ignorable and a 
general state n may be written as n = (ni, 77,2, 1 — ni— 77,2); 
if the total density is say 1, the the density of the third is 
simply 1 — ni — ri 2 . Following Zylka, we will demonstrate 
that K is star-likeP 



n\ 

FIG. 1. The set K for no = (0,2/7, 5/7), depicted in m-n 2 
space. The convex hull ch(K) is illustrated by the pink re¬ 
gion. Extreme points of ch(K) are indicated with red dots 
and labeled by the by the series of transformations Bij re¬ 
quired to reach them. Interior points are labeled with small 
black dots. K itself is the union of the red and black dots. 
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FIG. 2. The three lines covering K are depicted for no = 
(0,2/7, 5/7), emphasizing the star-shaped nature of K. 


To see this, first let us consider the example considered 
earlier, and previously^ namely the case of TV = 3 case 
with initial data no = (0, 2/7, 5/7) and e = (0,1,4). 

One can generate the entire set K by applying arbi¬ 
trarily long sequences of the transforms B^ to the ini¬ 
tial point no- We denote the set of points generated by 
k arbitrary Bjj the k th generation of K. The first six 
generations of K are plotted in Figure 1, overlaid with 
ch(K). Except for the original point no, the entirety of 
K lies along three distinct line segments joining pairs of 
extreme points, as depicted in Figure 2 . This arrange¬ 
ment is obvious if one considers that every state must 
have rii = nj, for some i and j, if the last transform 
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FIG. 3. K superimposed on contours of Wd = e • n (labeled 
by energy). The minimum Wd is found at the extreme point 
tiqB 22 ,B\ 2 ,B 12 — (3/8, 3/8,1/4) (labeled by a red circle in Figs. 
1 and 2) labeled here by a large white circle. 


applied in the sequence arriving at that state is Bij. 

The six extreme points (excepting no) can then be 
paired according to their final Bij. For example, the ex¬ 
treme points no5i2^i3 and noF> 23 ^i 3 are so paired, and 
the line joining them satisfies n\ = = 1 — n\ — n 2 , 

or ri 2 = 1 — 2n\. It is now clear that K must be 
star-shaped with respect to e = (1/3,1/3,1/3) because 
the three lines joining pairs of extreme points satisfy 
ni = n 2 = 1 — rii — ri 2 thereF 8 Likewise, ch(K ) for a 
three-level system is defined by at most seven extreme 
points (see Appendix A). Each generation approaches e 
more closely, and no extreme points are generated af¬ 
ter the third generation. (Consider the effect of Bij on 
a state n = (ni,ri 2 ,1 — n\ — n^). Then the Euclidean 
distance s of n from e satisfies s 2 = (n — e) • (n — e). 
Comparing the distances of n and nBij from e, one finds 
(n — e) • (n — e) — (nBij — e ) • ( nB^ — e) = \{ni — rij) 2 > 0. 
Thus each B\j brings a state closer to e.) 

Note that the set K is non-convexf^ since not all of 
the points in ch(K ) can be reached by sequences of B\j 
applied to no- Per the fundamental theorem of linear 
programming, the minimum Wd over the convex feasible 
region ch(K ) is found either at a single extreme point of 
ch(K) or along an edge of ch(K) joining two or more ex¬ 
treme points of ch(K). (In the latter case, the minimum 
Wd multiplicity is equal to the number of collinear ex¬ 
treme points.) However, the extreme points of ch(K) are 
all contained in K. Thus the linearity of the objective 
function Wd = e-n ensures that it has the same minimum 
over K and ch(K). 

Having identified the set of extreme points of ch(K ), 
£(ch(K )), in Fig. 1, the minimum Wd may be computed 
as min{e • £(ch(K))}. Fig. 3 overlays contours of the 
state energy e-n on iC, thereby illustrating the scheme 
and identifying the minimum energy state. In this case, 


the optimum level populations are nf = noB 2 3 BisBi 2 = 
(3/8, 3/8,1/4), yielding a minimum Wd = 11/8, as pre¬ 
viously calculated by exhaustive search. 

The sequence which minimizes the plasma energy can 
always be identified with this algorithm. Crucially, be¬ 
cause we have not assumed any ordering of the level ener¬ 
gies e, the slope of the state energy contours is arbitrary 
in the general case (cf. Fig. 3). Therefore, each of the 
seven extreme points identified is a possible solution to 
the energy minimization problem. 

IV. PARTIAL RELAXATION 

The discussion so far has assumed complete pairwise 
equalization of levels at every step, corresponding to 
transforms with a = 1/2. From a physical standpoint, 
densities could also be relaxed through a diffusion pro¬ 
cess that only partially equalizes population levels, corre¬ 
sponding to the case of T-transforms with 0 < a < 1/2. 

The a = 1/2 transforms B\j can be reached by repeat¬ 
edly applying a particular T-transform with 0 < a < 1 
formed from the same Qij , since for all 0 < a < 1, 
lim n _ >00 T’* = \(/ + Q^) = B^. Thus, the a = 1/2 
state space is in the closure of the state space for any 
fixed 0 < a < 1. At the same time, the 0 < a < 1/2 
state space is contained entirely in the convex hull of the 
a = 1/2 state space. To see this, simply observe that for 
0 < < 1/2, 

Tij = (1 — ol)I + OiQij = (1 — 2 ot)l + 2aBij, (2) 

where 0 < 2a < 1. Thus any 0 < a < 1/2 transform (on 
a state in the a = 1/2 state space) results in a state that 
is a convex combination of the initital state and the state 
resulting from the corresponding a = 1/2 transform, i.e. 
an interior point of ch(K). See Fig. 4 for a depiction of 
these possibilities. For any with 0 < a < 1/2, the 
new state is constrained to lie on the dashed line joining 
p and the corresponding pBij . 

It follows that it is sufficient to consider only a = 1/2 
transforms to identify the extreme points of the 0 < a < 
1/2 state space and accordingly solve the energy mini¬ 
mization problem. 

V. STRATEGIES FOR FINDING THE OPTIMUM 
SEQUENCE 

Two strategies for determining the sequence of trans¬ 
forms B^ required to attain the minimum -Wd state were 
conjectured. 2 Reproduced here, these strategies are: 

1. Diffusion of particles first between similar popula¬ 
tion levels, all other things being equal, eventually 
releases more energy. 

2. Depleting of particles the higher energy level first, 
all other things being equal, eventually releases 
more energy. 
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FIG. 4. Result of partial relaxations (T^ with 0 < a < 1/2) 
applied to a point p in the convex hull of K for no = 
(0,2/7, 5/7), as before. Partial relaxations result in interior 
points constrained to lie along the dashed lines joining p to 
each pBij. 


These strategies were only surmised based on trial- 
and-error experience. However, it turns out that they 
are intimately related to the more precisely put Propo¬ 
sition 2 , proposed and proven by Thon and Wallace! 9 It 
is worthwhile here to make use of the permutation of the 
set of initial populations, as introduced in Ref. |9] Given 
no = (^i, 712 , 723 ), the starting permutation is {1,2,3} 
if n\ < U 2 < 77 - 3 , as assumed throughout this work. 
Application of Bij exchanges the level numberings, e.g. 
if no has permutation {1,2,3}, noB \2 has permutation 
{2,1,3}. For a more detailed exposition of these per¬ 
mutations and emergent combinatorial techniques, please 
refer to Appendix B. 

Put briefly, Prop. 2 states that extreme points can be 
obtained by equalizing only adjacent level pairs. That 
is, any sequence of B\j resulting in an extreme point will 
effect only adjacent transpositions in the population per¬ 
mutation. More concretely, 770 H 12 and 770 H 23 are both 
extreme points, but no£> 13 , which averages the nonadja- 
cent first and third levels, is not. However, no£> 12 B \3 is 
an extreme point because the later application of Bis 
equalizes two adjacent levels, viz. no {1,2,3} 
noBi2 ~ {2,1,3} —noHi 2 Hi3 ~ {2,3,1}. 

Because the initial populations are assumed ordered 
increasing and the global minimum energy state will al¬ 
ways be located at an extreme point, any correct first 
step must be consistent with Strategy 1: the B\j chosen 
must mix two adjacent levels. 

It is clear that there are ‘dead ends’ among the pos¬ 
sible B{j sequences, where a state is reached with level 
populations decreasing with level energy, such that no 
more energy can be extracted. Any such stopping state 
has the level population permutation which is the re¬ 
verse of the energy level permutation. For example, given 


e = (ei,e 2 ,e 3 ) with e 2 < ei < e 3 ~ {2,1,3}, the stop¬ 
ping permutation is {3,1,2}, such that ns < n\ < 712 . 
n ^ {3,1, 2} can be reached from 77o r^j {1,2,3} by the 
sequence 770 H 23 H 13 . 

The three-level system has two extreme points with 
permutation {3,2,1}, the stopping permutation if the 
level densities and energies are both ordered increasing. 
In one case, the two lower-energy levels are first diffused; 
in the other, the two higher-energy levels are first dif¬ 
fused. The state energies of the extreme points cannot 
be ordered using a reduced set of variables (i.e. only some 
of the initial populations or energy levels), so the correct 
first step can only be determined in retrospect after a 
full calculation. The multiplicity of extreme points with 
a possible stopping permutation indicates that Strategy 
2 is not generally applicable. 

A further question of interest is the usefulness of an 
‘annealing’ strategy whereby a diffusion operation heats 
the system, resulting in a state with extra energy. Subse¬ 
quent B^ would lower the system energy to its minimum 
value, presumably smaller than the minimum value pos¬ 
sible without annealing. However, such a strategy cannot 
obtain a global minimum Wd . Note that any step which 
heats the system results in an inversion of the population 
permutation with respect to the correct stopping permu¬ 
tation. In order to obtain the stopping permutation, a 
subsequent diffusion operation is required on the two lev¬ 
els involved the heating step. 

As an example, consider a case mentioned previously 
with level energy permutation e ^ {2,1,3}, initial pop¬ 
ulation permutation Tio ^ {1,2,3}, and stopping permu¬ 
tation n rsj {3,1,2}. Note that the initial permutation 
contains two inversions with respect to the stopping per¬ 
mutation. Applying B 12 to the initial state results in the 
system absorbing energy, 

e • n 0 B 12 - e ■ n 0 = I(n 2 - ni)(ci - e 2 ). (3) 

By construction, each difference on the right hand side of 
Eq. ©> is positive and the system absorbs energy. Now 
consider the number of inversions: 770^12 results in the 
level population permutation n ~ {2,1,3}, containing a 
total of three inversions with respect to the stopping per¬ 
mutation {3,1,2}. The extra inversion changes the parity 
of the permutation and requires at least one additional 
diffusion operation to correct. In particular, the same dif¬ 
fusion operation B\ 2 must be repeated at some point in 
the relaxation process. Per Prop. 3 of Ref. 9] such a 
sequence does not generate an extreme point of ch(K ), 
even if both instances of the repeated B^ correspond to 
adjacent transpositions. However, there exists at least 
one extreme point with the stopping permutation. Be¬ 
cause K is star-shaped, any annealing strategy results 
in a non-optimal interior point (a convex combination of 
the uniform distribution e and an extreme point with the 
stopping permutation). Because the parity and number 
of inversions characterize any finite-length permutation, 
this conclusion holds in the general TV-level case. 
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VI. COMPLEXITY 

Given length-TV initial data, there are at most (^) 
possible states after one transformation has been ap¬ 
plied (fewer in case of degeneracy in the initial level 
populations). The second generation contributes up to 
N(N — 1) (TV — 2) (TV + 5)/8 unique states ? 1 Although the 
number of possible unique states grows rapidly with the 
passing generations, these later states spiral quickly to¬ 
ward the uniform distribution e = (1/TV, 1/TV,..., 1/TV). 
In fact, we can safely restrict our attention to the first 
(^) generations, a nevertheless enormous set for large TV. 

Thon and Wallacd^ presented an algorithm for gener¬ 
ating the extreme points of ch(K). In fact, the algorithm 
identifies the reduced (minimum length) sequences of ad¬ 
jacent transpositions (i.e. permutations exchanging only 
neighboring elements in the set) leading to each permu¬ 
tation in the TV th symmetric group Sn- For example, if 
TV = 3, the group £3 contains TV! = 6 possible permu¬ 
tations. The permutation {3,2,1} can be reached from 
the initial X = {1,2, 3} by two distinct sequences of adja¬ 
cent transpositions: (1, 2)(2, 3)(1, 2) and (2, 3)(1, 2)(2, 3). 
The other five possible permutations of X have unique 
reduced decompositions in adjacent transpositions; the 
total number of such decompositions for all of the permu¬ 
tations in £3 is therefore seven . 12 As noted, each possible 
permutation in Sn is a stopping state for an appropriate 
permutation of e. 

Any sequence of diffusion operations with a repeated 
Bij is not minimal and results in a permutation accessible 
with (^) or fewer B^. This follows because the reverse 
permutation, the longest permutation when expressed in 
adjacent transpositions, requires precisely (^) transposi¬ 
tions (corresponding to the number of inversions in the 
final permutation). For example, nB 2 3 B\ 3 B\ 2 B 2 3 has 
the permutation {2,3,1}, which could also be reached 
with nB^Bis. 

In the N = 3 case, the minimization of Wd is 
straightforward. There are Q) = 3 unique B^: 

F?i 3 5 andF? 23 - Without loss of generality, the ini¬ 
tial data no is increasing, so the initial permutation is 
{1, 2, 3}. The set K is covered by ch(K), which has seven 
extreme points: n 0 , n 0 B 12 , n 0 B 23 , n 0 B 12 B 13 , n o^ 23 ^i 3 , 
^ 0 ^ 12 ^ 13 ^ 23 , and noB 23 Bi 3 Bi 2 . 

There are e for which each of the seven extreme points 
can serve as the optimum system configuration. As 
noted, e ^ {1,2,3} is a special case in which the de¬ 
composition has maximum length (^) =3; there are 
two unique decompositions of the stopping permuta¬ 
tion {3,2,1} in adjacent transpositions. Equivalently, 
there are two extreme points with the stopping per¬ 
mutation: the minimum system energy is then Wd = 
min{e • noB 23 Bi 3 Bi 2: e • noBi 2 Bi 3 B 23 }. This is the most 
general scenario for the free energy optimization prob¬ 
lem, in which multiple extreme points with the correct 
stopping permutation must be compared. 

The problem has been reduced to evaluating the func¬ 
tion Wd = e • n at a finite number of known extreme 


points. Unfortunately, the upper bound on the num¬ 
ber of extreme points with stopping permutations is 
0(N N )PHH Because there is no a priori means of or¬ 
dering these points in energy, the exponential depth of 
the state tree is intrinsic to the optimization problem?^ 
Thus calculating the system energy accessible with dis¬ 
crete diffusive exchanges is NP-hard. 

VII. CONCLUSIONS 

By utilizing the (mostly ignored) literature on pairwise 
relaxation transformations, it is a relatively simple mat¬ 
ter to apply the techniques of convex optimization and 
combinatorics to answer the questions posed. In partic¬ 
ular, we find that for the three-level case, only seven se¬ 
quences need be considered. One further conclusion that 
can be drawn is that each of the sequences that must be 
considered will definitely be the solution at least for one 
energy vector. In general, a large but finite number of 
sequences must be checked to solve the energy extrac¬ 
tion problem. Moreover, partial relaxation or annealing 
strategies are never useful. 

The exposition here has also laid the foundation for 
analyzing the extremal properties of more complicated 
systems. One might make use of the results about K 
to extremize nonlinear objective functions, as would be 
necessary in e.g. systems with electrostatic self-energy. 
Alternately, formulating the optimization problem over 
a system of countable particles (as in integer program¬ 
ming) could shed light on the manipulation of degenerate 
matter. These systems and others like them could reveal 
the physical significance of the non-convexity of the state 
space. 

Although we have shown that identifying the maximal 
energy extraction solution is impractical in most cases 
of interest, the results obtained along the way nonethe¬ 
less provide useful constraints on diffusive schemes that 
should narrow the search for efficient strategies. 


ACKNOWLEDGMENTS 

The authors are very grateful for discussions with L. 
Shi and S. Davidovits. Particular thanks are due to 
Professor D. Tannor for critical discussions at an early 
stage of this work. Work supported by DOE Contract 
No. DE-AC02-09CH11466 and DOE NNSA SSAA Grant 
No. DE274-FG52-08NA28553. Oneofus(NJF) acknowl¬ 
edges the hospitality of the Weizmann Institute of Sci¬ 
ence, where he held a Weston Visiting Professorship dur¬ 
ing the time over which this work was initiated. 

1 N. J. Fisch and J.-M. Rax, Phys. Rev. Lett. 69, 612 (1992). 

2 N. J. Fisch and J.-M. Rax, Phys. Fluids B 5, 1754 (1993). 

3 N. J. Fisch, Fusion Sci. Technol. 65 (2014). 

4 M. C. Levy, S. C. Wilks, M. Tabak, S. B. Libby, and M. G. 

Baring, Nat. Commun. 5, 4149 (2014). 

5 C. S. Gardner, Phys. Fluids 6, 839 (1963). 




7 



FIG. 5. Geometric construction of K. Extreme points of 
ch(K) (yellow dots) depicted with lines rii = n, (solid) and 
lines joining extreme points (dashed). 


6 I. Y. Dodin and N. J. Fisch, Phys. Lett. A 341 , 187 (2005). 

7 F. Horn, in Proceedings of the 3rd European Symposium on 
Chemical Reaction Engineering (Pergamon, 1964) pp. 1-10. 

8 C. Zylka, Theor. Chim. Acta 68, 363 (1985). 

9 D. Thon and S. W. Wallace, Soc. Choice Welfare 22, 447 (2004). 

10 A. W. Marshall, I. Olkin, and B. C. Arnold, Inequalities: Theory 
of Majorization and Its Applications , 2nd ed. (Springer, New 
York, 2009). 

1:l OEIS Foundation Inc., “The On-Line Encyclopedia of Integer 
Sequences,” http://oeis.org/A239568 (2014), [Online; accessed 
13-Aug-2015]. 

12 OEIS Foundation Inc., “The On-Line Encyclopedia of Integer 
Sequences,” http://oeis.org/A246865 (2014), [Online; accessed 
13-Aug-2015]. 

13 OEIS Foundation Inc., “The On-Line Encyclopedia of Integer 
Sequences,” http://oeis.org/A005118 (2011), [Online; accessed 
13-Aug-2015]. 

14 R. P. Stanley, Europ. J. Combinatorics 5, 359 (1984). 

15 M. Yannakakis, J. Comput. Syst. Sci 43 , 441 (1991). 


Appendix A: Extreme point geometry 

Thon and Wallace proved 9 that their algorithm gener¬ 
ates the full set of extreme points for an arbitrary num¬ 
ber of energy levels N. However, the N = 3 case admits 
a straightforward geometric proof. All states resulting 
from the application of a Bij he along the line rii = rij. 
This constrains the possible paths through state space. 
In particular, the slope of the line joining any two points 
Pi and P2 such that P2 = PiBij is —1 if (i,j) = (1,2), 0 
if {hi) = (1,3), or oo if (i,j) = (2,3). 

Fig. [5] plots the extreme points and the three lines cov- 
ering K as well as these constrained trajectories. One 
extreme point is given (n, the initial data), and the other 
six are generated from choosing one of two routes through 
state space. Note that extreme points are only obtained 
by moving directly to a line rii = np without crossing 


any other such line. The lines are color-coded. Thus, for 
example, the partial relaxation along a red dotted line 
terminates on the red solid line; the red solid indicates 
complete relaxation of levels n\ and 712 , while all the 
states traversed along the dotted line are reachable. It 
can be seen from this geometrical representation that no 
extreme points are obtained by revisiting the same line: 
each diffusive transformation B^ reduces the distance of 
the state from the uniform distribution e. Thus, a state 
that revisits a line would necessarily be interior to the 
first state reached on that line. 


Appendix B: Combinatorial methods 

Thon and Wallac^ 9 simplified the structures of many 
proofs by making use of the combinatorial features of the 
Dalton transfer problem. In particular, they introduced 
the permutation of the individuals whose incomes were to 
be redistributed. In their notation, the initial permuta¬ 
tion is denoted v = (v(l),v(2 ),..., v(N)), where v(i) = k 
is interpreted to mean that the individual numbered k is 
the i th poorest. This representation is equivalent to the 
common two-line notation : 

f 1 2 ... N \ 

y v(l) v(2) • • • v(N) J 5 

where the first row lists indices of the permutation and 
the second row identifies the individuals. Note that the 
ordering of the columns is immaterial. The first column 
could be read as “individual v(l) is the poorest, listed 
first.” 

In their problem, as in ours, it is assumed that the in¬ 
dividuals (levels) are ordered initially from poorest (least 
populated) to richest (most populated), yielding an ini¬ 
tial permutation v(i) = i. Because the diffusion oper¬ 
ations treated in this work transfer level densities but 
leave the level energies unchanged, we lose no generality 
by specializing to such a permutation, which we denote 
{1,2,..., N }, reflecting nj < < • • • < Uq . (This is 

the ordered arrangement or one-line representation of the 
level density permutation.) In concordance with Ref. [9j 
we describe the levels nearest in population as neighbors , 
e.g. U 2 is neighbors with n\ and only. It turns out 
that the extreme points of ch(K ) can be reached only by 
sequences of B^ averaging such neighboring levels (Prop. 
2 of Ref. E|. 

Hoping for greater simplicity, we discuss sepa¬ 
rately the permutation of the level energies, denoted 
{w(l),w(2 ),... ,w(N)}, with w(i) = k meaning the level 
with the k th lowest energy has initially the i th smallest 
population. For example, if there is a complete popula¬ 
tion inversion, the highest-energy levels are initially most 
populated and e ^ {1,2,..., N}. If instead N = 3 and 
the second-highest energy level is most populated, fol¬ 
lowed in turn by the lowest and highest energy levels, 







one has e ~ {3,1, 2}, or, in two-line notation: 

(3 2 1 Wl 2 3\ 

V 2 1 3 7 V 3 1 2/ 

(Recall our assumption i = v(i), so that the first row 
of the two-line representation of the energy permutation 
identifies the states ordered by initial level populations. 
Thus the columns of the two-line representation could be 
parsed “greatest density in second-highest energy level, 
second-highest density in the lowest energy level, and 
least density in highest-energy level.”) 

Apart from equalizing the populations of levels i and 
j, one can consider the effect of the Bij on these permu¬ 
tations. As mentioned previously, the have no effect 
on the permutation of level energies, which remains fixed 
throughout the problem. However, each B\j changes the 
permutation of level populations, with the result that 
the numberings are exchanged, so that e.g. noBis has 
permutation {3, 2,1}. Considering the problem in n\-ri 2 
space, one realizes that the effect of a B\j is to move the 
system state to another of the N\ = 6 cells in ch(K ), 
each corresponding to a specific permutation of the level 
populations (cf. Fig. 1 of Ref. |9j) . The cells are sepa¬ 
rated by rays from e to six of extreme points in ch(K) 
(no excluded). 

One of the key insights of this approach is the connec¬ 
tion between population inversions and inversions in the 
permutation of level densities vis-a-vis the energy permu¬ 
tation. Using only the Bij : it is possible to reorder the en¬ 
tire set of level populations, so that the final distribution 


is decreasing with energy (we term the permutations cor¬ 
responding to such distributions stopping permutations). 

In fact, there is a unique stopping permutation in a 
given problem: the reverse of the energy permutation, 
i.e. {w(N),w(N — 1),..., w( 1)}, obtained from the com¬ 
position of w with the order-reversing permutation. Let 
w r denote the reverse of the permutation w. In two-line 
notation, one has 

1 2 ••• N \ f 1 2 ••• N\ 

w( 1) w( 2) ••• w(N) ) \N N-l ..• 1 ) 


f N 

N-l ■■ 

' 1 i 

1 2 •• 

■ N\ 

\w(N) 

w(N-l) •• 

• •Mi); 

\N N-l ■■ 

• lj 


( 1 2 • • • N \ 

y w(N) w(N — 1 ) ... w( 1 ) J ' 

Note that the reverse operation is an involution, such 
that (w r ) r = w. For example, {1,2,3} and {3,2,1} are 
reverse permutations. In our previous example, the level 
energy permutation was {3,1,2} such that 6% < e\ < 62 
for e = ( 61 , 62 , 63 ) and no = (n^n^ng). The reverse of 
the e permutation is {2,1,3}. Thus any stopping state 
has U2 < ni <713, such that the final population in the 
highest-energy level (‘level 2 ,’ when ordered by initial 
population) is smallest, and so on. 

By identifying the level energies’ permutation and its 
reverse, the search for the correct sequence of B^ can be 
greatly narrowed because only extreme points with the 
stopping permutation need be considered. 


